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ABSTRACT 

We use the large cosmological Millennium Simulation (MS) to construct the first all- 
sky maps of the lensing potential and the deflection angle, aiming at gravitational lens- 
ing of the CMB, with the goal of properly including small-scale non-linearities and 
non-Gaussianity. Exploiting the Born approximation, we implement a map-making 
procedure based on direct ray-tracing through the gravitational potential of the MS. 
We stack the simulation box in redshift shells up to z ^ 11, producing continuous 
all-sky maps with arcminute angular resolution. A randomization scheme avoids rep- 
etition of structures along the line of sight and structures larger than the MS box 
size are added to supply the missing contribution of large-scale (LS) structures to the 
lensing signal. The angular power spectra of the projected lensing potential and the 
deflection-angle modulus agree quite well with semi-analytic estimates on scales down 
to a few arcminutes, while we find a slight excess of power on small scales, which we 
interpret as being due to non-linear clustering in the MS. Our map- making proce- 
dure, combined with the LS adding technique, is ideally suited for studying lensing 
of CMB anisotropies, for analyzing cross-correlations with foreground structures, or 
other secondary CMB anisotropies such as the Rees-Sciama effect. 
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1 INTRODUCTION 

The cosmic microwave background (CMB) is character- 
ized both by primary anisotropies, imprinted at the last 
scattering surface, and by secondary anisotropies caused 
along the way to us by density inhomogeneities and re- 
scatterings on electrons that are freed during the epoch of 
reionization, and heated to high temperature when massive 
structures virialize. One of the interesting effects that can 
generate secondary anisotropies is the weak gravitational 
lensing of the CMB, which arises from the distortions in- 
duced in the geodesies of CMB photons by gradients in 
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the g r avitational matter potent ial (|Bartelmann fc Schneideil 
l200ll : iLewis fc Challinoj [iooil '). Forthcoming CMB probes 
do have the sensitivity and expected instrumental perfor- 
mance which may allow a detection of the lensing distor- 
tions of the primary CMB anisotropies, which would then 
also provide new insights and constraints on the expan- 
sion history of the univers e and on the process of cosmo- 
logical structur e formation (|Acauaviva fc Baccigalupill2006l : 
IHu et aEll2006l ). However, accurate predictions for the ex- 
pected anisotropies in total intensity and polarization are 
clearly needed for analyzing this future data, which demands 
for detailed simulated maps. 

The increasing availability of high-resolution N-body 
simulations in large periodic volumes makes it possible to di- 
rectly simulate the CMB distortions caused by weak lensing 
using realistic cosmological structure formation calculations. 
This work represents a first step in that direction. Existing 
studies already give access to statistical properties of the 
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Figure 1. Sketch of the adopted stacking and randomization process. The passage of CMB photons through the dark matter distribution 
of the Universe is followed by stacking the gravitational potential boxes of the MS, which are 500/i~^Mpc on a side (comoving). Shells 
of thickness 500/i~^Mpc are filled with periodic replicas of the box. All boxes (squares) that fall into the same shell are randomized with 
the same coordinate transformation (rotation and translation), which, in turn, differs from shell to shell. 



expected all-sky CMB lensing signal, such as the two-point 
correlation function and the power sp ectruni of th e lensing 
potential and deflection angle, see e.g. l|Lewisll2005l ) and ref- 
erences therein. This is based on 'semi-analytic' calculations 
that use approximate parameterizations of the non-linear 
evolution of the matter power spectrum. On the other hand, 
up to now N-body numerical simulations have been used to 
lens the CMB only on small patches of the sky in order to 
expl oit the practicality o f the flat-sky approximation, see 
e.g. (|Amblard et aZ.I [20041 ) and references therein. However, 
our approach of propagating rays through the forming dark 
matter structures gives access to the full statistics of the 
signal, including non-linear and non-Gaussian effects. Fur- 
thermore, it allows the accurate characterization of correla- 
tions of CMB lensing distortions with the cosmic large-scale 
structure, and with other foregrounds such as the Sunyaev- 
Zeldovich and Rees-Sciama effects. Hopefully this will allow 
improvements in the methods for separating the different 
contributions to CMB anisotropies in the data, which would 
be of tremendous help to uncover all the cosmological infor- 
mation in the forthcoming observations. 

From an experimental point of view, the improved preci- 
sion of the CMB observations, in particular that of the next 
generation experiment^, may in fact require an accurate de- 
lensing methodology and a detailed lensing reconstruction. 
CMB experiments targeting for instance the CMB polariza- 
tion, and in particular the curl component of the polariza- 
tion tensor, the so called B-modes from cosmological grav- 
itational waves, may greatly benefit from a precise knowl- 
edge of the lensing effects in order to separate them from 
the primordial cosmological signal (|Seliak fc Hiratal 12004 ). 
In particular, for a correct interpretation of the data from 
the forthcoming Planck satellite!, it will be absolutely es- 
sential to understand and model the CMB lensing, as the 
satellite has the sensitivity and overall instrumental perfor- 
mance for measuring the CMB lensing with good accuracy. 
We note that a first detection of CMB lensing in data from 
the Wilkinson Microwave Anisotropy Probe (WMAFjfl) com- 



^ See lambda.gsfc.nasa.gov for a complete list of operating and 
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bined with complem enta ry data has already been claimed by 
(jSmith et fflZ.ll200it ) and (Hirata et a;.ll2008l ). 

In this study we introduce a new methodology for the 
construction of all-sky lensing-potential and defiection-angle 
maps, based on a very large cosmo logical simulation, the 
Millennium run (|Springel et al\\2QQ± . As a first step in the 
analysis of the maps produced using the MS dark matter 
distribution, we have determined the interval of angular 
scales on which these maps match the semi-analytical ex- 
pectations, since we expected a lack of lensing power on 
large scales, due to the finite volume of the N-body simula- 
tion. To compensate for this effect, we have implemented a 
method for adding large-scale power which allows to recover 
the correct lensing signal on the scales outside this inter- 
val, i.e. on scales larger than the MS box size. At the other 
extreme, at the smallest resolved scales, we are interested 
in the question whether our maps show evidence for extra 
lensing power due to the accurate representation of higher- 
order non-linear effects in our simulation methodology. On 
these small scales, the impact of non-Gaussianities from the 
mapping of non-linear lenses is expected to be largest. 

This paper is organized as follows. In Section (5] we 
briefiy describe the basic aspects of lensing relevant to our 
work. In Section [31 we describe the N-body simulation and 
the details of our map-making procedure. In Section we 
present the lensing-potential and deflect ion- angle maps, and 
study the distribution of power in the angular domain. In 
Section [S] we provide a summary and discussion. 



2 LENSED MAPS OF THE CMB VIA THE 
BORN APPROXIMATION 

In what follows we will consider the small-angle scattering 
limit, i.e. the case where the change in the comoving sepa- 
ration of CMB light-rays, owing to the deflection caused by 
gravitational lensing from matter inhomogeneities, is small 
compared to the comoving separation of the undeflected 
rays. In this case it is sufficient to calculate all the rele- 
vant integrated quantities, i. e. the so-called lensing-potential 
and its angular gradient, the deflection-angle, along the un- 
deflected rays. This small-angle scattering limit corresponds 
to the so-called "Born approximation" . 

We treat the CMB last scattering as an instantaneous 
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process and neglect reionization. Adopting confo rmal time 
and comoving coordinates in a flat geometry (|Ma et al\ 
1 19951 ). the integral for the projected lensing-potential due 
to scalar perturbations with no anisotropic stress reads 



■ dr , 



*(n) = -2 



while the corresponding deflection-angle integral is 



a(n) 



r <i>(rn; 770 



■ dr , 



(1) 



(2) 



where r is the comoving distance, r« ~ 10'' Mpc is its 
value at the last-scattering surface, rjo is the present con- 
formal time, <1> is the physical peculiar gravitational po- 
tential generated by density perturbations, and [l/r]Vft is 
the two dimensional (2D) transverse derivative with respect 
to the liii e-of-sight p ointing in the direc t ion fi = ip) 
(IHuI I2OOOI: IFa rtclmani^ fc Schneideil I2OOII : iRefregierl l2003l : 
iLewis fc Chaninoill2006f r 

Actually, the lensing potential is formally divergent ow- 
ing to the 1/r term near r = 0; nonetheless, this divergence 
affects the lensing potential monopole only, which can be set 
to zero, since it does not contribute to the deflection-angle. 
In this way the remaining multipoles tak e a finite value and 
the le nsing potential field is well defined (|Lewis fc Challinoil 
I2OO6I ). Analytically, the full information about the deflection 
angle is contained in the lensing potential, but numerically 
the two equations ([T]) and @ are generally not equivalent, 
and it will typically be more accurate to solve the integral 
([2} directly to obtain the deflection angle instead of finite- 
differencing the lensing potential. 

If the gravitational potential "I> is Gaussian, the lensing 
potential is Gaussian as well. However, the lensed CMB 
is non-Gaussian, as it is a second order cosmological 
effect produced by cosmological perturbations onto CMB 
anisotropies, yielding a finite correlation between different 
scales and thus non-Gaussianity. This is expected to be 
most important on small scales, due to the non-linearity 
already present in the underlying properties of lenses. 

The most advanced approach develop ed so far fo r the 
construction of all-sky lensed CMB maps (Lewis 200i) em- 
ploys a se mi-analytical mod eling of the non-linear power 
spectrum l|Smith et al\ |2003| ). and derives from that the 
lensing potential and deflection angle templates assuming 
Gaussianity. This approach is therefore accurate for what 
concerns the two point correlation function of the lensing 
potential, as long as the non-linear two-point power of the 
matter is modeled correctly, but it ignores the influence of 
any statistics of higher order, which is expected to become 
relevant on small scales, where the non-linear power is most 
important. The use of N-body simulations to calculate the 
lensing has the advantage to possess a built-in capability of 
accurately taking into account all the effects of non-linear 
structure formation. On the other hand, the use of N-body 
simulations also faces limitations due to their limited mass 
and spatial resolution, and from their flnite volume, as we 
will discuss later on in more detail. 

For what concerns the line-of-sight integration in 
Eqs. Il} and ([2|, the Born-approximation along the unde- 
flected photon path holds to good accuracy and allows to 
obtain results which include the non-linear physics. Even on 



small scales, in fact, this approximation can be exploited in 
the small-angle scattering limit, i.e. for ty pical deflections 
being of the order of arcminutes or less (IHirata fc SeliakI 
I2OO3I : IShapiro fc Cooravll2006l ). For example, a single clus- 
ter typically gives deflection angles of a few arcminutes, 
while smaller structures, such as galaxies, lead to arcsecond 
deflections. Furthermore, it can be shown that the Born- 
approximation also holds in 'strong' lensing cases, provided 
that the deflection angles are equally small. Finally, second 
order corrections to the Born approximation (for instance 
a non- vanishing curl component) are expected to be sub- 
dominant with respect to the non-linear structur e evolution 
effects on small scales l|Lewis fc ChallinoJlioO^ . For these 
reasons, we argue that this approximation should be accu- 
rate enough for calculating all-sky weak lensing maps of the 
CMB based on cosmological N-body simulations. 



3 MAP-MAKING PROCEDURE FOR THE 
MILLENNIUM SIMULATION 

The Millennium Simulation (MS) is a high-resolution N- 
body simulation ca rried out by the Virgo Consortium 
(jSpringel et al\\200^ ). It uses N = 2160^ ~ 1.0078 X 10^" 
collisionless particles, with a mass of 8.6 X lO*/i"'M0, to 
follow structure formation from redshift z — 127 to the 
present, in a cubic region 500 /i~' Mpc on a side, and with 
periodic boundary conditions. Here h is the Hubble con- 
stant in units of 100kms~'Mpc~'. With ten times as many 

f iarticles as the previous largest computations of this kind 
Colberg ef aZll2000l : lEvrard et aZll2002l : IWambsganss et al\ 



), it features a substantially improved spatial and time 
resolution within a large cosmological volume. 

The cosmological parameters of the MS are as follows. 
The ratio between the total matter density and the criti- 
cal one is = 0.25, of which ilt = 0.045 is in baryons, 
while the density of cold dark matter (CDM) is given by 
ficDM = — f^i). The spatial curvature is assumed to be 
zero, with the remaining cosmological energy density made 
up by a cosmological constant, Ha ~ 0.75. The Hubble 
constant is taken to be Ho — 73kms~'Mpc~'. The pri- 
mordial power spectrum of density fluctuations in Fourier 
space is assumed to be a simple scale-invariant power law of 
wavenumber, with spectral index = 1. Its normalization 
is set by the rms fluctuations in spheres of radius 8h~^ Mpc, 
(jg = 0.9, in the linearly extrapolated density field at the 
present epoch. The adopted parameter values are consistent 
with a combined analysis of the 2dF Galaxy Redshift Sur- 
vey (2dfGRS) and the fir st year WMAP data (|Colless et al\ 
l200ll : ISpergel et aZ.llioO^ . 

Thanks to its large dynamic range, the MS has been 
able to determine the non-linear matter power spectrum 
over a larger range of scales than possible in earlier works 
(jjenkins et al\ Il998l l . Alm ost five orders of m agnitude in 
wavenumber are covered ijSpringel et al\\200^ . This is a 
very important feature for studies of CMB lensing, as we ex- 
pect that this dynamic range, combined with the method for 
adding large-scale structures described in the next section, 
can be leveraged to obtain access to the full non-Gaussian 
statistics of the lensing signal, limited only by the the max- 
imum angular resolution resulting from the gravitational 
softening length and particle number of MS. We stress again 
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Figure 2. The simulated all-slcy map of the lensing potential computed with the map-malting procedure combined with the LS adding 
method as described in the text. 



that the lensed CMB is non-Gaussian even if the underly- 
ing lenses do possess a Gaussian distribution. Moreover, the 
non-linear evolution of large scale structures produces a de- 
gree of non-Gaussianity in the lenses distribution which con- 
tributes to the non-Gaussian statistics of the lensed CMB 
on small scales. This non-Gaussian contribution can be com- 
puted only via the use of N-body simulations which are able 
to accurately describe the non-linear evolution of the lenses. 
These non-linearities are known to alter the lensed tempera- 
ture power spectrum of CMB anisotropics by about ~ 0.2% 
at ^ ~ 2000 and by ~ 1% or more on smaller scales. But, 
much more notably, they introduce ~ 10% co rrections to 
the B-mode polarization power on all the scales (|Lewisll2005l : 
iLewis fc Challinoill2006l ). 

Our map-making procedure is based on ray-tracing of 
the CMB photons in the Born approximation through the 
three-dimensional (3D) field of the peculiar gravitational po- 
tential. The latter is precomputed and stored for each of 
the MS output times on a Cartesian grid with a mesh of 
dimension 2560^ that covers the comoving simulation box 
of volume (500 /i~^Mpc)''. The gravitational potential itself 
has been calculated by first assigning the particles to the 
mesh with the clouds-in-cells mass assignment scheme. The 
resulting density field has then been Fourier transformed, 
multiplied with the Green's function of the Poisson equa- 
tion in Fourier space, and then transformed back to obtain 
the potential. Also, a slight Gaussian smoothing on a scale 
Ts equal to 1.25 times the mesh size has been applied in 
Fourier space in order to eliminate residual anisotropies on 
the scale of the mesh, and a deconvolution to filter out the 
clouds-in-cells mass assignment kernel has been applied as 
well. The final potential field hence corresponds to the den- 
sity field of the MS (which contains structures down to the 



gravitational softening length of 5/i ^kpc) smoothed on a 
scale of ~ 200/i~^kpc. 

In order to produce mock maps that cover the past light- 
cone over the full sky, we stack the peculiar gravitational po- 
tential grids around the observer (which is located at 2 = 0), 
producing a volume which is large enough to carry out the 
integration over all redshifts relevant for CMB lensing. For 
simplicity, we only integrate out to z, = 11.22 in this study, 
which corresponds to a comoving distance of approximately 
r, ~ 7236/i~^Mpc with the present choice of cosmological 
parameters. Indeed, the lensing power from still higher red- 
shifts than this epoch is negligible for CMB lensing, as we 
will discuss in the next section. But we note that our method 
could in principle be extended to still higher redshifts, up to 
the starting redshift 2 — 127 of the simulation. 

The above implies that the simulation volume needs to 
be repeated roughly 14.5 times along both the positive and 
negative directions of the three principal Cartesian axes x,y, 
and 2, with the origin at the observer. However, the spacing 
of the time outputs of the MS simulation is such that it cor- 
responds to an average distance of 140/i~^Mpc (comoving) 
on the past light-cone. We fully exploit this time resolution 
and use 53 outputs of the simulation along our integration 
paths. In practice this means that the data corresponding 
to a particular output time is utilized in a spherical shell of 
average thickness 140/i~^Mpc around the observer. 

The need to repeat the simulation volume due to its fi- 
nite size immediately means that, without augmenting large- 
scale structures, the maps will suffer from a deficit of lensing 
power on large angular scales, due to the finite MS box size. 
More importantly, a scheme is required to avoid the repe- 
tition of the same structures along the line of sight. Previ- 
ous studies that constructed simulated light-cone maps for 
small patches of the sky typically simply randomized each 
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Figure 3. Simulated all-sky maps of the defloction-angle components along the i9 direction (top panel), and along the ip direction (bottom 
panel), in radians. 



of the repeated boxes along the past lightcone by apply- 
ing independent ran dom translations and reflections (e.g. 
ISpringel et all 1200 j ). However, in the present application 
this procedure would produce artefacts like ripples in the 
simulated deflect ion- angle field, because the gravitational 
field would become discontinuous at box boundaries, leading 
to jumps in the defiection angle. It is therefore mandatory 
that the simulated lensing potential of our all sky maps is 
everywhere continuous on the sky, which requires that the 
3D tessellation of the peculiar gravitational potential is con- 
tinuous transverse to every line of sight. 



Our solution is to divide up the volume out to 2* into 
spherical shells, each of thickness 500 /i^^Mpc comoving (ob- 
viously the innermost shell is actually a sphere of comoving 
radius 250 /i~^Mpc, centered at the observer). All the simu- 
lation boxes falling into the same shell are made to undergo 
the same, coherent randomization process, i.e. they are all 
translated and rotated with the same random vectors gener- 
ating a homogeneous coordinate transformation throughout 
the shell. But this randomization changes from shell to shell. 
Figure [1] shows a schematic sketch of this stacking process. 
For simplicity, the diagram does not illustrate the additional 
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Figure 4. Top panel: Simulated all-slcy map of the deficction-anglc modulus (in radians), obtained with the map-making procedure 
combined with the LS adding method as described in the text. 



shell structure stemming from the different output times 
of the simulation. As discussed before, this simply means 
that the underlying potential grid is updated on average 
3-4 times with a different simulation output when integrat- 
ing through one of the rotated and translated 500/i~^Mpc 
shells, but without changing the coordinate transformation. 
Notice that our stacking procedure eliminates any preferred 
direction in the simulated all-sky maps. 

In order to define the gravitational potential at each 
point along a ray in direction n, we employ spatial tri- 
linear interpolation in the gravitational potential grid. It is 
then easy to numerically calculate the integral potential for 
each ray, based for example on a simple trapezoidal formula, 
which we use in this study. Obtaining the deflection angle 
could in principle be done by finite differencing a calculated 
lensing potential map, either in real space or the harmonic 
domain. However, the accuracy of this approach would de- 
pend critically on the angular resolution of the map. Also, 
the sampling of the gravitational potential in the direction 
transverse to the line-of-sight varies greatly with the dis- 
tance from the observer, so in order to extract the maximum 
information from the simulation data down to the smallest 
resolved scales in the potential field, we prefer to directly 
integrate up the defiection angle vector along each light ray 
in our map. For this purpose we first use a fourth-order 
finite-differencing scheme to compute the local 3D grid of the 
gradient of the gravitational potential, which is then again 
tri-linearly interpolated to each integration point along a 
line-of-sight. In this way, we calculate the deflection angle 
directly via equation ([5} along the paths of undeflected light 
rays. 

Finally, we need to select a pixelization of the sky with a 
set of directions n = {■&, ip). We here follow the standard ap- 



proach introduced by the HEALPi>Q hierarchical tessellation 
of the unit sphere (|Gorski et a/.ii2005. ). 



4 SIMULATED MAPS OF THE LENSING 
POTENTIAL AND DEFLECTION ANGLE 

In Figs. [51 [31 and we show full-sky maps of the lens- 
ing potential, the deflection angle i9/(^-components, and the 
deflection angle modulus \cx\, respectively, obtained with 
the map- making technique described in the previous section 
combined with a semi-analytic procedure (to be explained 
below) augmenting the lensing power on scales beyond the 
MS box size. These maps are generated with a HEALPix pix- 
elization paramet er A^aide = 2048, an d have an angular reso- 
lution of ~ 1.72' iGorski et ttZ.ll2005l ). with 50331648 pixels 
in total. 

Several interesting features should be noted in these 
maps. The distribution of the lensing potential, where the 
monopole and dipole have been cut to simplify the visual in- 
spection, appears to be dominated by large features, which 
are probably simply arising from the projection of the largest 
scale gravitational potential fluctuations along the line-of- 
sight. However, the strength of local lensing distortions in 
the CMB cannot be directly inferred from the map of the 
lensing potential, as for the lensing deflection only the gra- 
dient of the potential is what really matters. 

The maps showing the lensing deflection angle compo- 
nents have interesting features as well. First of all, the sig- 
nal in the two components of the deflection angle appears to 
possess two morphologically distinct regimes, characterized 
on one hand by a diffuse background distribution, caused 

healpix.jpl.nasa.gov 
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Figure 5. Top panel: The power spectrum of the simulated lensing potential map of Fig. |2] (blue solid line), compared with the power 
spectrum of the lensing potential obtained with the CAMB code (dashed-dotted black line), which also includes an estimate of the 
non-linear contributions ijSmith et adl20o3 ). The red dashed-Sdotted and orange dotted lines differ only in the starting redshift for the 
line-of-sight integration used in the map-making. While the result shown in red begins at z = 0, the orange line gives the result for a start 
at z = 0.22. Finally, the light-green dashed line represents the linear lensing potential power spectrum from the CAMB code. Bottom 
panel: The power spectrum (in radians squared) of the simulated deflection angle modulus map shown in the upper panel of Fig. 141 (blue 
solid line), compared with the power spectrum (dashed-dotted black line) of the synthetic deflection angle modulus map obtained as a 
Gaussian realization from the CAMB lensing potential power spectrum. As above, the red dashed-Sdotted and orange dotted lines differ 
only in the starting redshift of the line-of-sight integration, as labelled. The red line is for the full redshift interval, the orange one for a 
start at 2 = 0.22, as described in the text. 



probably by the lines-of-sight where no dominant structures 
are encountered, and on the other hand by sharp features, 
caused probably by massive CDM structures which give rise 
to the largest deflections in the line-of-sight integration it- 
self. The same features are evident in the map of the mod- 
ulus of the deflection-angle. 

The mean value of |q:| in our simulated maps is 2.36', 
while its standard deviation is 1.25'. The latter has to be 
compared with the corresponding value obtained via the an- 



gular differentiation of synthetic Gaussian maps produced 
with the lensing potential power spectrum generated by the 
publicly available Code for Anisotropies in the Microwave 
Background (CAMl|f|) using the MS cosmological parame- 
ters, as we explain in detail below. We find only a 0.03% 
difference for the rms of the |a|-maps from MS and CAMB, 



^ See camb.info. 
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when using the same maximum redshift of line-of-sight in- 
tegration, i.e. Zmax ~ 11.22. On the other hand, if we set 
Zmax ~ 1100 in CAMB, we find that our estimate is ~ 1.7% 
smaller than the semi-analytic one, due to the missed contri- 
bution from sources beyond z ~ 11 in our map- making pro- 
cedure. For comparison, we also evaluate semi-analytically 
the expected change in the standard deviation of |q:| when 
inserting in CAMB more recent estima tes of the cosmo- 
logical parameters ijKomatsu et al\\200^ . In this case the 
rms from MS is ~ 6% and ~ 4.2% greater than the semi- 
analytical prediction when in CAMB we set Zmax ~ 11.22 
and 2max ~ 1100, respectively. 

The lensing potential and deflection angle maps of 
Figs.[5]and|3]have been obtained combining the map- making 
procedure described in the previous section with the method 
for adding large-scale power that we now explain. 
Firstly, we have measured the power spectra of the simu- 
lated maps obtained from the MS scales only, i.e., using the 
routine ANAFAST of the HEALPix package, we have inde- 
pendently measured the power spectra of the lensing po- 
tential (C**) and deflection angle modulus (C"") of the 
MS simulated maps, without exploiting the relations be- 
tween the lensing potential and the ±l-spin components of 
the deflection angle, which hold in the spherical harmonic 
domain tHu,.2000. V Secondly, using the MS cosmological pa- 
rameters, we have evaluated the semi-analytical power spec- 
trum of the lensing potential from CAMB, incl uding the es- 
timate of the contribution from non-linearity (Smith et al\ 
120031 ') and stopping the line-of-sight integration redshift up 
to 2 = 11.22. Using the lensing potential power spectrum 
from CAMB, we have then produced the corresponding syn- 
thetic map (and its angular differentiation) obtained as a 
Gaussian realization generated with the HEALPix code SYN- 
FAST, in order to produce the synthetic map of the deflec- 
tion angle modulus from the semi-analytic expectations of 
CAMB. From this map we have then extracted the power 
spectrum of the deflection angle modulus, and after decon- 
volution from the HEALPix pixel window function, we have 
compared it, together with the lensing potential power spec- 
trum, to the corresponding deconvolved MS power spectra. 

The top panel of Fig. [S] shows the primary result of 
this comparison. The black dashed-dotted line represents 
the semi-analytic prediction of the lensing potential angular 
power spectrum obtained from CAMB as discussed above. 
This has been compared with the red dashed-3dotted line 
obtained with the map-making procedure previously de- 
scribed, and which represents the result for the full inte- 
gration starting at z = and ending at z — 11.22. In this 
case, a power deficiency on large scales with respect to the 
semi-analytical prediction is evident, and confined to a mul- 
tipole range corresponding to one degree or more in the sky. 
The same for the orange dotted line which gives the MS lens- 
ing potential power spectrum obtained from a line-of-sight 
integration starting at a redshift of z — 0.22 and ending 
at z — 11.22; comparing the two curves, a power decrease 
at low £ is easily observable in the orange dotted line, with 
respect to the red dashed-Sdotted one, illustrating the influ- 
ence of the lack of comoving scales greater than 500 /i~^Mpc 
in the MS. As expected, this effect is evident in the multi- 
pole range corresponding to a few degrees or more, which 
is about the size of the MS box at the redshift most rele- 
vant for CMB lensing, i.e. z ~ 1. However, towards larger 



£, the deficit of large-scale power quickly decreases, and be- 
comes negligible at scales I > 350. Between these scales and 
I ~ 2500, there is quite good agreement between the MS 
lensing power spectrum and the semi-analytic prediction, 
but at 2500 < I < 4000 the full MS signal for the lensing 
potential actually slightly exceeds the semi-analytic result. 
On this multipole range, the red dashed line is dominated 
by Poisson noise, but the slight excess of power is clearly 
observable from the orange dot line, in which there is no 
contribution from the low redshift integration at z < 0.2. 
We ascribe this power excess to the matter non-linearities 
accurately reproduced from the Millennium Simulation. Fi- 
nally, at I ~ 4000 the MS signal is dominated by Poisson 
sampling noise from low-redshift potential integration. In 
fact, at very low redshifts, the 1.72' angular resolution of 
our map is comparable and even smaller than the intrin- 
sic angular resolution corresponding to the spatial grid of 
the 3D gravitational potential field we use. This is evident 
in Fig. [S] where we compare the map angular resolution of 
1.72' (red dashed line) with the effective angular resolution 
corresponding to the intrinsic grid spacing (195/i~^kpc) of 
the 3D gravitational potential field as function of redshift. 
Because the line-of-sight integral for the projected lensing 
potential involves a 1/r weighting term, the resulting noise 
terms are unfavourably amplified when the lensing potential 
is considered. 

The comparison above has been used to evaluate the 
multipole range, < Z < 350, not covered by the MS scales. 
On this interval we have applied the LS adding method: 
from the CAMB and MS maps of the lensing potential, we 
have extracted the two corresponding ensembles ^t;^*^^ 
and ^i^^ of spherical harmonic coefficients, respectively. 
Since on low multipoles the effects of the non-Gaussianity 
from the non-linear scales are negligible and the 'i'lm are 
independent, we have generated a joined ensemble of "^im, 
where = for < « < 350 and = 'ffi^' for 

I > 350. Finally, we have generated the synthetic maps of 
the lensing potential and deflection angle as non-Gaussian 
constrained realizations, inserting the as input in SYN- 
FAST, as shown in Figs. [^Bl 

These maps have the peculiarity of reproducing the 
non-linear and non-Gaussian effects of the MS non-linear 
dark matter distribution at multipoles I > 350, while at the 
same time including the contribution from the large scales 
at I < 350, where the lensing potential follows mostly the 
linear trend as shown from the light-green dot-dashed line 
in Fig. [S] The blue solid curve in the same Figure represents 
the resulting power spectrum of the lensing potential map 
after the LS addition. 

The bottom panel of Fig. [S] shows the corresponding 
power spectra for the physically and numerically more 
meaningful deflection angle. Here we show a comparison of 
the power spectrum of the deflection angle modulus {C"") 
measured for the MS simulated maps, in the absence of LS 
supplying, with the semi-analytic prediction constructed 
with CAMB and SYNFAST, as explained above. Again, we 
flnd a deficit of power on large scales, and a reassuring 
agreement over about one order of magnitude in I on 
intermediate scales. However, a slight excess of power over 
the semi-analytic predictions is easily seen at I > 2500. 
As previously mentioned, it can be attributed to the 
non-linear evolution of the MS structures. Finally, the blue 
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solid line represents the power spectrum extracted from 
the deflection angle modulus map of Fig. after adding 
large-scale structures. 

Our map making procedure offers very good resolution 
at the most important redshift for lensing of the CMB, 
~ 1 (see also Fig. [71 , where the intrinsic angular resolution 
of our potential grid is six times better than the angular 
resolution of the full-sky map. We therefore think that 
this higher small-scale power is a direct result of the more 
accurate representation of non-linear structure formation 
in our map simulation methodology. In fact, in our current 
maps we are still far from probing the most non-linear 
scales accessible in principle with our simulation. Those are 
a factor 40 smaller (namely 5/i~^kpc) than resolved by the 
potential grid we have employed. However, using such a fine 
mesh is currently impractical, and would lead to angular 
resolutions in full-sky maps that are unaccessible even by 
the Planck satellite. However, for a smaller solid-angle 
of the map, these sc ales can be probed with a different 
ray-tracing technique (jHilbert et aZ.|[2007l ). 




Figure 6. Comparison between the angular resolution of 1.72' of 
our full-sky maps (red dashed line) and the redshift-dependent, 
effective angular resolution (blue solid line) corresponding to the 
intrinsic grid spacing (~ 200 ?i~^kpc) of the three-dimension grav- 
itational potential field constructed from the Millennium Simula- 
tion. 



We note that the semi-analytic prediction for the power 
spectrum of the deflection angle modulus has been evalu- 
ated as an angular gradient in the harmonic domain of a 
synthetic lensing potential Gaussian map; that is accurate 
since in this approach we work with Fourier modes right 
from the start anyway. From a numerical point of view, the 
integral and derivative operators in Eq. ^ do however not 
commute, even if they analytically do, in the sense that fi- 
nite differencing our measured projected potential will not 
necessarily give the same result as numerically integrating 
the deflection angle along each line of sight. The latter ap- 
proach is more accurate, expecially at very high resolution, 
and it has been used by us in the comparison above since 
numerically integrating the deflection angle along each line 
of sight allows to preserve the contribution from the non- 
linear scales in a more efficient way than simply operating 
in the harmonic domain. 

Finally, we consider the distribution of the deflection 
angle power along the line-of-sight. In Fig. [T] we show the 
cumulative and differential variance of the deflection angle as 
a function of redshift. We see that the most important con- 
tributions to the final signal stem from z ^ 1, i.e. about half 
ways between the last scattering surface and the observer, 
as expected. This also allows us to assess the relative error 
introduced by stopping the integration at z ~ 11, which is 
of the order of a few percent, as mentioned above. 



5 CONCLUSIONS 

We constructed the first all sky maps of the cosmic mi- 
crowave background (CMB) weak-lensing potential and de- 
flection angle based on a high-resolution cosmological N- 
body simulation, the Millennium Run Simulation (MS). The 
lensing potential and deflection angle are evaluated in the 
Born approximation by directly ray-tracing through a three- 
dimensional, high-resolution mesh of the evolving peculiar 
gravitational potential and its gradient. The time evolution 
is approximated by 53 simulation outputs between redshift 



2 = and z ~ 11, each used to cover a thin redshift interval 
corresponding to a shell in the past light-cone around the 
observer. To prevent artificial repetition of structures along 
the line-of-sight, while at the same time avoiding disconti- 
nuities in the force transverse to a line-of-sight, we tessellate 
shells of comoving thickness corresponding to the size of the 
box (500ft~^Mpc) with periodic replicas which are coher- 
ently rotated and translated within each shell by a random 
amount. Moreover, in order to include the contribution to 
the lensing signal from the scales larger than the MS box 
size, we have implemented a method for adding large-scale 
structure as described in the text. 

Using the Hierarchical Equal Area Latitude Pixelization 
(HEALPIX) package for obtaining a uniform sky-coverage, 
we have constructed simulated CMB lensing maps with ~ 5 
million pixels and an angular resolution of ^ 1.72', based on 
potential fields calculated on 2560"^ meshes from the Millen- 
nium simulation. In the present study, we analyze the power 
spectrum of the lensing potential and the deflection an- 
gle, and compare it with predictions made by semi-analytic 
approaches. We note that our general approach for map- 
making can be extended to other CMB foregrounds, includ- 
ing the Integrated Sachs- Wolfe (ISW) and Rees-Sciama ef- 
fects at low redshifts, as well as estimates of the Sunyaev 
Zel'dovich (SZ) effects, or of the X-ray background. This 
will in particular allow studies of the cross-correlation of the 
lensing of CMB temperature and polarization with these ef- 
fects, which will be the subject of a forthcoming study. In 
our approach we do not take into account the contributions 
of the baryonic physics to the lensing effects on the CMB. 
We expect in fact that these contributions could be non- 
negligible only on the typical scales of cluster cores and be- 
low, thus well above I ~ 3000. Our comparison of the angular 
power spectrum of the lensing-potential and the deflection- 
angle with semi-analytic expectations reveals two different 
regimes in our results. First, for multipoles up to Z ~ 2500, 
our simulated maps produce a lensing signal that matches 
the semi-analytic expectation. Second, we flnd evidence for 
a slight excess of power in our simulated maps on scales cor- 
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Figure 7. The cumulative and differential variance of the deflec- 
tion angle map as a function of redshift. The symbols mark the 
different output times of the Millennium Simulation. The vertical 
dashed line gives the redshift that corresponds to the 50% quar- 
tile of the total variance in our maps, which is approximately at 
^ ~ 1. The dotted lines mark the 5% and 95% percentiles, indicat- 
ing that 90% of our signal in the deflection angle power spectrum 
is produced in the redshift range z 0.1 to ^ ~ 6.0. Note how- 
ever that we have lost a few percent of the total power due to our 
truncation of the integration at z = 11.22. When included, this 
would slightly shift these precentiles to higher redshift. 



responding to few arcminutes and less, which we attribute 
to the accurate inclusion of non-linear power in the Millen- 
nium simulation. It will be especially interesting to study 
the non-Gaussianities in the signal we found and its implied 
consequences for CMB observations. 

The new method proposed here demonstrates that an 
all-sky mapping of CMB lensing can be obtained based on 
modern high-resolution N-body simulations. This opens the 
way towards a full and accurate characterization of CMB 
lensing statistics, which is unaccessible beyond the power 
spectrum with the existing semi-analytical techniques. This 
is relevant in view of the forthcoming CMB probes, both as 
a way to detect, extract and study the CMB lensing signal, 
which carries hints on the early structure formation as well 
as the onset of cosmic acceleration, and as a tool to distin- 
guish CMB lensing from the Gaussian contribution due to 
primordial gravitational fluctuations. 
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